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摘 要 : 为 了 估计 给 水 系统 节点 流量 和 节点 水 压 参数 ， 实 时 模拟 给 水 系统 的 运行 状态 ， 本 文 建立 了 参数 
估计 最 小 二 乘 数学 模型 。 该 模型 是 带 等 式 、 非 负 约束 的 非 线性 优化 问题 。 利 用 Wolfe 既 约 梯度 法 
和 Rosen 投影 梯度 法 求解 优化 模型 ， 讨 论 了 模型 求解 步骤 ， 得 出 了 Rosen 投影 梯度 法 中 投影 矩 
阵 的 快速 计算 方法 。 编 制 了 求解 模型 的 计算 机 程序 ， 应 用 于 实际 给 水 系统 的 参数 估计 ， 估 计 的 节 
点 水 压 与 实测 节点 水 压 的 平均 误差 为 1.35m， 表 明 模 型 及 其 求解 方法 是 可 行 的 ， 且 Rosen 法 的 速 


度 更 快 。 
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1 ”给 水 系统 及 其 参数 估计 问题 


给 水 系统 是 由 管 段 和 节点 两 类 基本 元 素 组 成 的 水 流 输送 系统 ， 由 管 段 围 成 的 封闭 环 称 为 
基 环 ， 如 图 1 所 示 。 节 点 参数 包括 节点 流量 Q@ 和 节点 水 压 态 。 管 段 参数 包括 直径 D、 长 度 工 、 
粗粮 度 C、 流 量 g、 流 速 v 和 水 头 损 失 h， 其 中 D、LKL、C 是 静态 参数 ， 反 映 了 管 段 的 输 水 能 
力 ， 可 折算 成 一 个 参数 5， 称 为 摩 阻 。 管 段 流速 。、 水 头 损失 及 和 节点 水 压 矿 都 可 以 由 管 段 流 
景 g 和 摩 阻 S 直 接 求 出 ， 因 而 不 是 独立 的 参数 。 独 立 参 数 只 有 Q@、gq 和 5S。 各 参数 之 间 的 关系 
满足 以 下 几 个 方程 : 


96 (9) 27 (©) 
图 1: 给 水 系统 简 图 
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1) 管 段 方程 
管 段 的 水 头 损失 由 管 段 流量 和 摩 阻 计算 得 到 


j = 9797 ， j= 1,2,.…,P, (1) 


其 中 忆 为 管 段 总 数 ，nn 为 常数 。 
2) 节点 连续 方程 
流 进 节点 的 流量 等 于 流出 节点 的 流量 


》 aijgj = 一 Qi i=1,2,.…,N, (2) 
j=1 


其 中 入 为 节点 总 数 ，aij 为 系数 ， 其 取 值 方法 是 ， 当 管 段 了 与 节点 ?相连 且 流 入 节点 ;时 取 -1， 
流出 节点 i 时 取 1， 当 管 段 7 与 节点 i 不 相连 时 取 0。Qi 为 节点 流量 ， 流 出 为 正 ， 流 进 为 负 。 
公式 (2) 是 由 入 个 方程 组 成 的 方程 组 ， 写 成 矩阵 形式 为 


[4] : [a] = ~[Q]), (3) 


其 中 [4 为 由 系数 @ij 组 成 的 NN x 书 维 矩阵 ，[g] = [ga ,9qpl7，[Q] = [QQ ,Qwj 

方程 组 (2)、(3) 中 ， 任 何 一 个 方程 可 由 其 余 N 一 1 个 方程 相 加 求 得 ， 因 此 方程 组 (2)、(3) 中 
只 有 NN 一 1 个 独立 方程 。 

3) 环 能 量 方 程 

在 每 个 基 环 中 ， 顺 时 针 方 向 流动 的 管 段 水 头 损 失 之 和 ， 等 于 逆 时 针 方 向 流动 的 管 段 水 头 损 
失 之 和 ， 若 设 顺 时 针 方向 管 段 水 头 损 失 为 正 ， 逆 时 针 方 向 管 段 水 头 损失 为 负 ， 则 基 环 内 所 有 管 
段 水 头 损 失 代数 和 为 零 。 


P 
DLi(S;g7) =0, 1=1,2,..,L,, (4) 
1 


其 中 工 , 为 基 环 总 数 ，Lij 为 系数 ， 其 取 值 方法 是 ， 当 管 段 7 在 基 环 1 中 且 顺 时 针 方向 流动 时 
取 1， 道 时 针 方 向 流动 时 取 -1， 当 管 段 ?不 在 基 环 ! 中 时 取 0。 
上 式 是 由 工 。 个 方程 组 成 的 方程 组 ， 写 成 矩阵 形式 为 


四 [so =0， (5) 


其 中 [ 且 为 由 系数 Lj 组 成 的 Ls x 忆 维 矩阵 ，[Sg"] = [Si1g?, S2g2,:… ,Spg%]T。 

给 水 系统 的 节点 数 N、 管 段 数 P 和 基 环 数 工 存在 如 下 关系 : 已 = N++ 工 一 1， 则 方 
程 (3) 和 (5) 的 独立 方程 总 数 为 P 个 。 

若 给 水 系统 输入 的 流量 、 水 压 一 定 ， 并 且 已 知 其 余 节 点 的 输出 流量 时 ， 可 联 立 方 
程 (3)、(5) 式 ， 求 得 书 个 管 段 流 量 q9。 再 求 各 管 段 流 速 、 水 头 损失 和 每 个 节点 的 水 压 。 这 
一 计算 过 程 称 为 管 网 平 差 。 

水 压 、 水 量 是 保证 供水 服务 质量 的 基础 。 一 个 运行 良好 的 给 水 系统 应 能 够 为 用 户 提供 充足 
的 水 量 和 适宜 的 水 压 ， 同 时 使 供水 成 本 最 小 。 及 时 了 解 给 水 系统 各 节点 水 压 、 各 管 段 的 流速 、 
流量 、 水 力 坡 降 等 参数 是 非常 必要 的 ， 这 将 有 利于 指导 给 水 系统 的 优化 调度 ， 增 强 运行 管理 人 
员 的 信心 ， 发 现 管 网 中 存在 的 问题 ， 为 管 网 的 改 扩建 提供 依据 。 

由 于 给 水 系统 是 一 个 大 型 复杂 的 系统 ， 各 节点 输出 的 流量 @; 随 用 户 用 水 量 的 变化 而 变化 ， 
表现 出 很 强 的 时 变性 和 不 确定 性 ， 因 而 无 法 利用 (3)、(5) 式 计算 管 段 流量 及 节点 水 压 。 目 前 只 
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能 采用 很 粗略 的 方法 了 解 给 水 系统 的 运行 状态 ， 在 给 水 系统 中 安装 监测 仪表 ， 测 定 管 网 中 部 分 
节点 水 压 和 管 段 流量 ， 并 将 监测 数据 实时 传送 到 控制 中 心 。 由 于 监测 设备 投资 大 ， 而 管 网 节点 
数 众 多 ， 为 了 节省 投资 ， 只 能 实测 少量 的 数据 ， 因 而 对 系统 状态 的 了 解 很 不 全 面 。 

本 文 研究 的 参数 估计 方法 ， 利 用 实测 的 给 水 系统 部 分 管 段 流量 和 节点 水 压 数 据 ， 结 合 静态 
的 管 网 图 形 、 管 段 参数 数据 ， 反 推 估计 所 有 节点 的 用 水 量 、 节 点 水 压 和 管 段 流量 ， 达 到 实时 模 
拟 给 水 系统 运行 状态 的 目的 。 


2 ”给 水 系统 参数 估计 模型 的 建立 


给 水 系统 参数 估计 的 前 提 是 ， 已 知 管 网 的 图 形 结构 和 管 段 摩 阻 $， 已 知 某 时 刻 各 水 源 输入 
水 量 ， 实 测 了 管 网 中 部 分 节点 水 压 Ho 和 管 段 流量 qoj 。 参 数 估计 的 基本 思路 和 方法 是 ， 假 设 
各 节点 流量 @i， 使 其 总 和 等 于 各 水 源 供水 量 之 和 ， 利 用 假设 的 节点 流量 通过 (3)、(5) 式 平 差 
计算 得 所 有 管 段 流量 9%j， 进 而 求 得 所 有 节点 水 压 瑟 。 比 较 实 测 点 处 的 计算 值 及 、g; 与 实测 
值 Ho:、go; 的 差别 。 如 果 误 差 较 大 ， 说 明 假 设 的 各 节点 流量 与 实际 发 生 的 节点 流量 相差 较 大 ， 
应 调整 。 利 用 调整 后 的 节点 流量 再 进行 平 差 计 算 ， 直 至 计算 值 与 实测 值 比较 接近 ， 这 时 的 节点 
流量 接近 真实 的 值 ， 此 时 的 节点 流量 即 为 参数 估计 的 结果 ， 计 算 结 果 即 反映 了 给 水 系统 运行 状 
态 的 全 貌 。 

为 了 避免 在 进行 节点 流量 调整 时 的 盲 自 性 ， 减 小 节点 流量 搜索 范围 ， 缩 短 搜索 时 间 ， 对 节 
点 流量 的 搜索 做 一 定 的 约束 ， 使 节点 流量 围绕 某 一 个 值 搜索 。 理 由 是 ， 虽 和 然 节点 流量 是 随时 变 
化 的 ， 但 相对 而 言 ， 某 些 节点 的 流量 较 大 ， 而 另 一 些 节 点 流量 较 小 ， 相 互 间 存在 一 定 的 比例 关 
系 ， 这 种 关系 在 一 段 时 间 内 变化 不 是 很 大 ， 可 利用 以 往 的 用 户 水 表 记 录 的 用 水 量 数据 计算 各 节 
点 流量 的 比例 关系 。 每 次 进行 状态 模拟 时 ， 利 用 实测 的 总 供水 量 、 各 节点 流量 比例 关系 ， 粗 略 
计算 各 节点 流量 Q0;:， 将 粗略 计算 的 节点 流量 作为 节点 流量 的 搜索 中 心 。 由 于 粗略 计算 的 节点 
流量 可 能 存在 较 大 的 误差 ， 因 此 ， 只 能 让 它 对 节点 流量 搜索 产生 很 小 的 约束 力 。 

由 以 上 分 析 ， 建 立 如 下 给 水 系统 参数 估计 的 最 优化 模型 


NH Na N 
minf(Q) = > Wi(D)(Ho: — Hi)? + > Wo(j)(g0; — q;)* + > Wali)(Qoi ~ Qi)? (6) 
j=1 i=1 


l=1 


N 
st 站 9 0) 
Qi > 0, i=1,2,.…,N, 


其 中 Wi() 为 节点 测 压 参数 权重 ; W2( 让 为 管 段 测 流 参数 权重 ， Was() 为 节点 流量 约束 权 
重 ， NH 为 测 压 节点 数 ，N 为 测 流 管 段 数 ; NN 为 管 网 节点 总 数 ，Q 为 节点 流量 集合 ，Q = 
[Q1, Q2,… , QN]7; Qa 为 系统 总 供水 量 。 

测 压 点 节点 水 压 权 重 W 是 用 以 反映 测量 参数 的 准确 度 、 精 确 度 的 权重 系数 。 实 际 取 值 时 ， 
选 定 一 个 准确 度 、 精 确 度 最 高 的 测量 参数 ， 令 其 权重 为 1， 其 余 参 数 权重 根据 实际 情况 取 小 于 
等 于 1 的 值 。 节 点 流量 约束 权重 Ws(i) 应 取 较 小 的 值 ， 其 值 可 取 0.005 左右 ， 这 样 使 节点 流量 
粗略 计算 值 对 节点 流量 搜索 的 约束 很 小 。 


3 ”优化 模型 的 求解 
上 节 所 建立 的 优化 模型 (6)-(7) 是 一 个 具有 线性 等 式 约束 和 非 负 约束 的 非 线性 规划 问题 ， 可 


822 工 程 数 学 学 报 第 26 卷 


用 Wolfe 既 约 梯度 法 、Rosen 投影 梯度 法 或 罚 函数 法 求解 ， 因 前 两 者 针对 线性 约束 的 计算 效率 
较 高 ， 本 文采 用 前 两 者 。 

3.1 Wolfe 既 约 梯度 法 

3.1.1 求解 步骤 

该 法 首先 选 定 一 个 节点 流量 Qw 作为 基本 变量 四， 将 其 表示 成 其 余 入 -1 工 个 节点 流量 的 函数 


六 一 1 


Qw=Qu- > 9， 
?一 1 
从 而 消去 了 等 式 约束 条 件 。 
目标 函数 关于 变 景 Q; 的 既 约 梯度 为 
Es a Ee A 1,2,.……,N—1. 


OQ: OQN 6Qi OQ: OQN 
为 了 满足 非 负 约束 条 件 ， 对 搜索 方向 修正 如 下 


jw -neo ni(QW) <0, 
| -rn(Qw)， 当 ri(QW) > 0 


由 于 
人 一 1 
QN = Qa— > Qi 
t=1 
则 有 


于 是 得 到 节点 流量 迭代 公式 为 
i 过 QW) + tdit®), i= 2 wl; (8) 


其 中 心 (人 为 搜索 方向 ，t 为 搜索 步 长 ,上 > 0。 
为 了 使 新 搜索 的 节点 流量 Qt+D) 满足 非 负 约束 ，t 应 在 0 ~ tmax 之 闻 进 行 一 维 寻 优 ，timax 按 
下 式 计算 


(Kk) < o}. (9) 


ba min{ ~ 


Qt ) 
d; 

假设 各 节点 流明 ， 构 建 一 组 节点 流量 的 可 行 解 Q(， 使 其 总 和 等 于 总 供水 量 ， 进 行 半 差 计 
算 ， 计算 目标 函数 值 1。 若 f 值 较 大 ， 按 (8) 式 进 行 节点 流量 调整 。 

设 定 一 个 很 小 的 数 e。 > 0， 当 ao < = 时 ， 优 化 计算 终止 ， 此 时 的 解 集 @ 即 为 参数 估计 的 最 
优 解 。 用 节点 流量 @ 代 入 方程 组 (3)、(5)， 求 解 所 有 管 段 流量 g， 计 算 节点 水 奈 互 ， 得 到 给 水 
系统 运行 状态 模拟 结果 。 

3.1,2 梯度 计算 

对 i 求 一 阶 偶 导 数 
总 


荔 - -二 Oo 个 


->》 2W2(j)(go; — 中 ) 玛 A [ea ) 0h WE ee 
j=1 
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其 中 下、gj 扩 分 别 为 利用 Q() 平 差 计算 得 到 的 测 压 点 水 压 和 测 流 管 段 流量 ，Q,( 是 加 
的 一 个 分 量 。 

在 上 式 中 ， 比 较 难 计算 的 是 8H1/6Q; 和 Bq;/8Q;， 可 利用 (3)、(5) 式 计算 ,计算 方 法 参见 
文献 [2]。 

3.2 ”Rosen 投 影 梯度 法 

3.2.1 求解 步骤 

将 约束 条 件 (7) 写成 矩阵 形式 加 


[Dl: (Q] = Qa, (10) 

[B]: [Q] > 0， (11) 
其 中 [D] 3 [1, 1: , ;1]， 共 人 个 元 素 ; [Q] 一 [Q1, Q2,…… ,QnN]T; [B] 为 N x NN 维 单位 方 阵 。 
由 矩阵 [Dj 和 [B] 构 建 窍 阵 [C]， 将 [DJ 作为 [C] 的 第 一 列 。 假 设 在 优化 计算 过 程 中 ， 

在 (11) 式 表示 的 入 个 不 等 式 约束 中 ， 有 m 个 式 子 满足 等 式 条 件 ， 将 [B] 中 与 等 式 约束 对 应 


的 行 向 量 转 置 ， 作 为 [C] 的 另外 m 个 列 向 量 。 则 [OC] 为 N x (m 十 1) 维和 矩阵 。 
第 kk 十 1 步 节点 流量 的 迭代 公式 为 Q%1D) = Qt +tdto ，dto 的 最 优 搜索 方向 为 


d®) = —[PJ[Vf(Q™)], (12) 
其 中 do = [a 的 ,gd 罗 ),.…. ,dkOIT， [局 为 投影 矩阵 ， 是 入 阶 方 阵 ，[P] = I 一 [Cl([C]7[C])-1[O]7; 
of of of 
Viton = [og 30%' "505} 


搜索 步 长 上 同样 应 满足 (9) 式 的 条 件 。 

3.2.2 ”投影 矩阵 的 快速 计算 

投影 矩阵 [P] 计算 的 关键 是 矩阵 [C] 的 确定 ， 而 确定 [C] 的 关键 是 确定 其 后 mm 列 数据 。 理 论 
上 ，[Q] 的 后 m 列 数据 是 从 [B] 中 选 出 的 对 应 于 Qt = 0 的 那些 行 的 数据 。 但 是 在 实际 计算 过 
程 中 ， 随 着 迭代 的 进行 ， 某 些 节点 流量 Q@(5 只 是 趋 近 于 0， 但 并 不 可 能 等 于 0。 当 @Q(5 足够 小 
时 ， 由 (9) 式 计算 的 如 ax 将 很 小 ， 这 限制 了 和 迭代 计算 前 进 的 步伐 ， 使 欠 代 原 地 踏步 。 因 此 ， 
必须 设 定 一 个 较 小 的 值 Quin， 当 @ 人 9 小 于 Qmin 时 认为 它 近似 等 于 0， 并 将 该 节点 流量 对 应 
的 [Bj 中 的 那 行 元 素 设 为 [C] 中 的 一 列 元 素 。 经 过 这 种 近似 简化 后 ， 和 迭代 得 以 进行 ， 收 敛 速度 
较 快 。 在 后 续 的 计算 实例 中 采用 的 Qmin = 2L/s， 约 为 平均 节点 流量 的 3%。 

在 迭代 计算 过 程 中 ，{[C] 在 不 断 变化 ， 因 而 需要 不 断 地 计算 投影 矩阵 [ 媚 。 每 次 计算 投影 
矩阵 时 ， 需 要 计算 ([Cj7[C]) 的 逆 矩 阵 ， 计 算 量 很 大 。 经 过 理论 分 析 发 现 ，([CJ7IC]) 的 逆 矩 阵 
是 mm 十 1 维 方 阵 ， 可 直接 采用 下 式 ， 不 需 进 行 计 算 ， 使 计算 速度 大 幅度 提高 。 


6 -6 -6 .… -6 
-8 1+6 8 .1 

(CYTICN) =| -8 6 1+6 : 8 |, 
二 . a i 8 


-6B 8 Bb BB 1+8 


N—m—1 
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4 应 用 实例 


我 们 利用 上 述 方法 编制 了 给 水 系统 参数 估计 的 计算 机 程序 ， bd nde 
估计 。 某 大 城市 给 水 管 网 共有 167 个 节点 、277 个 管 段 、9 个 供水 水 源 ， 真 径 150mm 以 上 管 
累计 长 度 1220km， 供 水 面积 142km2。 对 各 段 管道 的 直径 、 长 度 、 粗 糙 度 等 结构 2 
查 测 定 ， 建 立 了 管 网 拓扑 结构 模型 。 为 了 了 解 给 水 系统 的 运行 状况 ， 在 某 一 时 刻 同时 测定 了 管 
网 中 115 个 节点 的 水 压 ， 并 测定 了 该 时 刻 各 水 源 供水 量 ， 总 供水 量 为 12076 L/s。 

按 文献 图 中 的 方法 ， 从 115 个 已 知 节点 水 压 中 选择 了 24 个 节点 水 压 ， 作 为 目标 函数 中 的 实 
测 节点 水 压 Ho:， 利 用 参数 估计 计算 程序 模拟 出 167 个 节点 的 水 压 。 将 模拟 节点 水 压 与 已 知 节 
点 水 压 进 行 了 对 比 ， 结 果 如 表 1。 


表 1: 给 水 系统 参数 估计 结果 (Rosen 法 ) {单位 : m) 


模拟 节 | 实测 节 | 水 压 模拟 节 | 实测 节 | 水 压 | 节点 | 模拟 节 
点 水 压 | 点 水 压 | 误差 点 水 压 | 点 水 压 | 误差 编号 | 点 水 压 


实测 节 
点 水 压 


1 440.75 = 站 
2 441.93 | 441.6 0.33 
3 441.77 | 439.5 2.27 
4 441.93 | 439.0 2.93 
5 442.55 
6 445.41 445.2 0.21 
7 444 444.0 0 

8 437.99 | 439.0 | -1.01 
9 438.74 | 437.1 1.64 
10 438.96 | 438.2 0.76 
11 439.14 | 439.3 | -0.16 
12 440.02 | 440.5 | -0.48 
13 440.32 | 440.2 0.12 
14 443.65 441.39 | 442.8 | -1.41 
15 436.22 442.98 | 442.4 0.58 
16 436.22 443.14 | 442.2 0.94 
17 442.08 443.47 = ee 
18 441.46 443.51 442.1 1.41 
19 437.7 443.55 > = 
20 442.65 444.18 | 445.2 | -1.02 
21 439.45 439.89 | 439.6 0.29 
22 439.14 441.6 440.7 0.9 
23 441.44 446.78 ss -一 
24 436.22 446.67 | 445.4 1.27 
25 436.94 445.53 一 
26 442.01 444.27 | 443.7 0.57 
27 440.93 438.19 | 437.6 0.59 
28 441.52 441.01 439.0 2.01 
29 439.24 全 一 

439.44 | 438.8 
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( 续 ) 表 1: 给 水 系统 参数 估计 结果 (Rosen 法 ) (单位 : m) 


点 | 模拟 节 | 实测 节 | 水 压 ‖ 节点 | 模拟 节 | 实测 节 于 节点 | 模拟 节 
号 | 点 水 压 | 点 水 压 | 误差 | 编号 | 点 水 压 | 点 水 压 | 误 点 水 压 
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两 种 方法 计算 结果 统计 见 表 2。 总 体 上 看 ， 水 压 误差 较 小 ， 能 满足 管 网 实时 模拟 的 要 求 。 两 


种 方法 的 计算 过 程 稳定 且 收 敛 ， 而 Rosen 法 计算 速度 更 快 。 


表 2: 参数 估计 计算 结果 统计 


平均 水 压 误差 (m) | 最 大 水 压 误差 (m) 二 
1.38 3.26 
1.35 2.96 


计算 机 CPU 速度 为 1.7G， 内 存 为 512M。 


计算 耗 时 *(s) 
96 


Wolfe 法 
Rosen 法 


5 ”结论 
本 文 研究 给 水 系统 的 参数 估计 问题 ， 主 要 结果 如 下 : 
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1) 建立 了 给 水 系统 参数 估计 的 优化 模型 ， 在 管 网 拓扑 结构 已 知 的 基础 上 ， 利 用 少量 的 实 
测 节点 水 压 和 管 段 流量 ， 估 计 所 有 节点 流量 ， 模 拟 计算 管 网 所 有 节点 的 水 压 ， 以 获得 给 水 系统 
运行 状态 的 全 貌 。 

2) 讨论 了 采用 Wolfe 既 约 梯度 法 和 Rosen 投影 梯度 法 求解 目标 函数 的 方法 ， 编 制 了 计算 程 
序 。 计 算 表 明 ， 两 种 方法 都 是 可 行 的 ， 而 Rosen 法 的 运行 速度 更 快 。 

3) 将 模拟 计算 方法 应 用 于 某 城市 给 水 系统 参数 估计 计算 ， 结 果 表 明 ， 估 计 的 精度 较 高 ， 
计算 过 程 稳定 收敛 ， 计 算 速 度 较 快 ， 能 满足 给 水 系统 参数 估计 的 需要 ， 为 给 水 系统 运行 状态 实 
时 模拟 提供 了 新 的 途径 。 
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Abstract: In order to simulate the operating state of a water supply system, a least-square mathematic 
model has been built to estimate node water flow and water pressure head in the water supply system. 
The model is a nonlinear optimization model with equality and no negative constraints. The Wolfe 
gradient and the Rosen projection algorithm are used to solve this model. The framework for the 
model solving has been discussed. A quick method to calculate the Rosen projection matrix has been 
determined. A program has been developed to solve the model, and is applied to estimate parameters 
of an actual water supply system. The average difference of node water head between the estimated 
values and the measured values is 1.35m. This shows that the model and its solving algorithm are 
feasible. Moreover, we find that the Rosen projection algorithm is quicker than the Wolfe gradient 
algorithm. 

Keywords: water supply system; parameter estimate; optimization model; Wolfe gradient; Rosen 
projection 


